#include "Matrix.h"

using namespace std;

Matrix::Matrix(){
	Matrix(0,0,CPU_MEM);
}
Matrix::Matrix(string file_name){
	ifstream file;
	
	file.open(file_name.c_str(), ios::in | ios::binary);
	if(!file){
		cerr<<"Unable to open file: "<<file_name<<endl;
		exit(2);
	}
	file.read((char*)(this),sizeof(*this));
	memory=CPU_MEM;
	data = new REALN[ndata];
	file.read((char*)data,size);
	file.close();
}
Matrix::Matrix(UINTN r,UINTN c,short mem){
	row=r;
	col=c;
	ndata=row*col;
	esize=sizeof(REALN);
	size=ndata*sizeof(REALN);
	if(mem==CPU_MEM){
		data=new REALN[ndata];
	}else{
		CUDA_SC(cudaMalloc((void**)&data,size));
	}
	mode = 'n';
	memory=mem;
	ext_memory=false;
}
Matrix::Matrix(UINTN r,UINTN c,short mem,short mmode){
	row=r;
	col=c;
	ndata=row*col;
	esize=sizeof(REALN);
	size=ndata*sizeof(REALN);
	if(mem==CPU_MEM){
		data=new REALN[ndata];
	}else{
		CUDA_SC(cudaMalloc((void**)&data,size));
	}
	mode = mmode;
	memory=mem;
	ext_memory=false;
}
Matrix::Matrix(UINTN r,UINTN c,REALN* new_data,short mem){
	row=r;
	col=c;
	ndata=row*col;
	esize=sizeof(REALN);
	size=ndata*sizeof(REALN);
	data=new_data;
	//__asm__("int3");
	/*
	if(mem==CPU_MEM){
		data=new_data;
	}else{
		CUDA_SC(cudaMalloc((void**)&data,size));
		CUDA_SC(cudaMemcpy(data,new_data,size,cudaMemcpyHostToDevice));
	}*/
	mode = 'n';
	memory=mem;
	ext_memory=true;
}
Matrix::Matrix(Matrix* m,short mem){
	row=m->getRow();
	col=m->getCol();
	ndata=row*col;
	esize=sizeof(REALN);
	size=ndata*sizeof(REALN);
	memory=mem;
	mode = m->getMode();
	
	if(mem==CPU_MEM && m->getMemorySpace()==CPU_MEM){
		data=new REALN[ndata];
		memcpy((void *)data,(void *)m->getData(),size);
	}else if(mem==CPU_MEM && m->getMemorySpace()==DEVICE_MEM){
		data=new REALN[ndata];
		CUDA_SC(cudaMemcpy(data,m->getData(),size,cudaMemcpyDeviceToHost));
		//__asm__("int3");
	}else if(mem==DEVICE_MEM && m->getMemorySpace()==CPU_MEM){
		CUDA_SC(cudaMalloc((void**)&data,size));
		CUDA_SC(cudaMemcpy(data,m->getData(),size,cudaMemcpyHostToDevice));
	}else if(mem==DEVICE_MEM && m->getMemorySpace()==DEVICE_MEM){
		CUDA_SC(cudaMalloc((void**)&data,size));
		CUDA_SC(cudaMemcpy(data,m->getData(),size,cudaMemcpyDeviceToDevice));
	}
	ext_memory=false;
}
Matrix::~Matrix(){
	//cout<<ext_memory<<endl;
	if(ext_memory){
		return;
	}
	if(memory==CPU_MEM){
		delete(data);
	}else{
		//MatrixMath::ShowMem();
		CUDA_SC(cudaFree(data));
		//MatrixMath::ShowMem();
		//cout<<endl;
	}
}
void Matrix::save(string file_name){
	ofstream file;
	REALN* data_temp;
	
	file.open(file_name.c_str(), ios::out | ios::binary);
	if(!file){
		cerr<<"Unable to open file: "<<file_name<<endl;
		exit(2);
	}
	file.write((char*)(this),sizeof(*this));
	if(memory==CPU_MEM){
		file.write((char*)data,size);
	}else{
		data_temp = new REALN[ndata];
		CUDA_SC(cudaMemcpy(data_temp,data,size,cudaMemcpyDeviceToHost));
		file.write((char*)data_temp,size);
		delete(data_temp);
	}
	file.close();
}
	
	
void Matrix::setUni(REALN v){
	UINTN i,j;
	Matrix* temp;
	if(memory==CPU_MEM){
		temp=this;
	}else{
		temp = new Matrix(row,col,CPU_MEM,mode);
	}
	for(i=0;i<row;i++){
		for(j=0;j<col;j++){
			temp->set(v,i,j);
		}
	}
	if(memory==DEVICE_MEM){
		CUDA_SC(cudaMemcpy(data,temp->getData(),size,cudaMemcpyHostToDevice));
		delete(temp);
	}
}
void Matrix::setGrow(REALN v,REALN s){
	UINTN i,j;
	REALN t;
	Matrix* temp;
	
	if(memory==CPU_MEM){
		temp=this;
	}else{
		temp = new Matrix(row,col,CPU_MEM,mode);
	}
	t=v;
	for(i=0;i<row;i++){
		for(j=0;j<col;j++){
			temp->set(t,i,j);
			t+=s;
		}
	}
	if(memory==DEVICE_MEM){
		CUDA_SC(cudaMemcpy(data,temp->getData(),size,cudaMemcpyHostToDevice));
		delete(temp);
	}
}
UINTN Matrix::getRow(){
	return row;
}
UINTN Matrix::getCol(){
	return col;
}
UINTN Matrix::getNData(){
	return ndata;
}
UINTN Matrix::getESize(){
	return esize;
}
UINTN Matrix::getSize(){
	return size;
}
REALN* Matrix::getData(){
	return data;
}
char Matrix::getMode(){
	return mode;
}
short Matrix::getMemorySpace(){
	return memory;
}
void Matrix::transpose(){
	UINTN t;
	t=row;
	row=col;
	col=t;
	if(mode=='n'){
		mode='t';
	}else{
		mode='n';
	}
	/*if(!(row==1||col==1)){
		cout<<"transposition is allowed only for vectors"<<endl;
		exit(2);
	}*/ //+++
}
void Matrix::swapdim(){
	UINTN t;
	t=row;
	row=col;
	col=t;
}
void Matrix::setMode(char m){
	mode=m;
}
UINTN Matrix::getLeading(){
	if(mode=='n'){
		return col;
	}else{
		return row;
	}
}
bool Matrix::operator!=(Matrix& other){
	UINTN i,j;
	
	if(col!=other.getCol() || row!=other.getRow()){
		return true;
	}
	for(i=0;i<row;i++){
		for(j=0;j<col;j++){
			if(get(i,j)!=other.get(i,j)){
					return true;
			}
		}
	}
	return false;
}
REALN Matrix::operator%(Matrix& other){
	UINTN i,j;
	REALN max=0.0;
	REALN error;
	
	if(col!=other.getCol() || row!=other.getRow()){
		return true;
	}
	for(i=0;i<row;i++){
		for(j=0;j<col;j++){
			if(get(i,j)!=other.get(i,j)){
				if((get(i,j))==0){
					error =fabs((double)other.get(i,j));
				}else{ 
					error=fabs(((double)(get(i,j)-(double)other.get(i,j)))/
								((double)(get(i,j))));
				}
				if(error>1){
					cout<<error<< "\t\t"<<get(i,j)<<"\t\t"<<other.get(i,j)<<endl;
					max+=1.0;
				}
				/*
				if(error>max){
					max=error;
					cout<<max<<" "<<i<<" "<<j<<" "<<get(i,j)<<" "<<other.get(i,j)<<endl;
				}
				*/
			}
		}
	}
	return max;
}
bool Matrix::operator==(Matrix& other){
	return !(*this!=other);
}

void Matrix::print(){
	UINTN i,j;
	Matrix* temp;
	
	if(memory==CPU_MEM){
		temp=this;
	}else{
		temp = new Matrix(this,CPU_MEM);
	}
	
	for(i=0;i<row;i++){
		for(j=0;j<col;j++){
			cout<<temp->get(i,j)<<"\t";
		}
		cout<<endl;
	}
	cout<<endl;
	
	if(memory==DEVICE_MEM){
		delete(temp);
	}
}
void Matrix::print(UINTN m){
	UINTN i;
	Matrix* temp;
	
	if(memory==CPU_MEM){
		temp=this;
	}else{
		temp = new Matrix(this,CPU_MEM);
	}
	for(i=0;i<m;i++){
		cout<<temp->getData()[i]<<endl;
	}
	if(memory==DEVICE_MEM){
		delete(temp);
	}
}
void Matrix::printInfo(){
	cout<<row<<"x"<<col<<" - elements: "<<ndata<<", total size: "<<
				size<<", mode: "<<mode<<", memory: ";
	if(memory==CPU_MEM){
		cout<<"CPU";
	}else{
		cout<<"GPU";
	}
	cout<<endl;
}
void Matrix::printInfo(string name){
	cout<<name<<" = ";
	printInfo();
}
void Matrix::append(Matrix* add,UINTN start,UINTN end){
	UINTN i;
	
	for(i=0;i<add->getRow();i++){
		if(add->getMemorySpace()==CPU_MEM){
			//cout<<i<<endl;
			//cout<<data<<" "<<data+size<<" "<<data+i*col+start<<endl;
			//cout<<add->getData()<<" "<<add->getData()+add->getSize()<<" "<<add->getData()+i*add->getCol()+end<<endl;
			/*for(j=0;j<end-start;j++){
				data[0]=add->getData()[i*add->getCol()+start+j];
			}*/
			memcpy(&(data[i*col+start]),&(add->getData()[i*add->getCol()]),(end-start)*sizeof(REALN));
		}else{
			CUDA_SC(cudaMemcpy(&(data[i*col+start]),&(add->getData()[i*add->getCol()]),
						(end-start)*sizeof(REALN),cudaMemcpyDeviceToHost));
		}
	}
}

	/*
	if(mode==ROW_WISE){
		for(i=0;i<ndata;i++){
			if(i%col==0 && i!=0){
				cout << endl;
			}
			cout << data[i] <<"\t";
		}
		cout << endl <<endl;
	}else{r
		for(i=0;i<col;i++){
			for(j=0;j<row;j++){
				cout << data[i+j*col] <<"\t";
			}
			cout << endl;
		}
		cout << endl <<endl;
	}*/
	
void Matrix::printRaw(UINTN l){
	UINTN i;
	//cout<<mode<<endl;

	for(i=0;i<l;i++){
		cout<<data[i]<<endl;
	}
	cout<<endl<<endl;
}
void Matrix::fast_compare(string fn1,string fn2,REALN max_error,REALN round_to_zero){
	ifstream f1,f2;
	UINTN i,j;
	UINTN t=0;
	Matrix* m1;
	Matrix* m2;
	REALN error;
	REALN min;
	
	min = 1/(pow(10,round_to_zero));
	
	m1=new Matrix(fn1);
	m2=new Matrix(fn2);
	
	for(i=0;i<m1->getRow();i++){
		for(j=0;j<m1->getCol();j++){
			if(m1->get(i,j)<min && m2->get(i,j)<min){
				continue;
			}
			if(m2->get(i,j)==0){
				error = fabs(m1->get(i,j));
			}else{ 
				error = fabs((m1->get(i,j)-m2->get(i,j))/m2->get(i,j));
			}
			if(error>max_error){
				cout<<i<<" "<<j<<" "<<"\t"<<m1->get(i,j)<<" "<<m2->get(i,j)<<endl;
				t++;
			}
		}
	}
	cout<<t<<endl;
	f1.close();
	f2.close();
}


